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Recently, cluster methods have been used to solve a variety of sign problems including those that arise in the 
presence of fermions. In all cases an analytic partial re-summation over a class of configurations in the path 
integral was necessary. Here the new ideas are illustrated using the example of QCD at a finite density of static 
quarks. In this limit the sign problem simplifies since the fermionic part decouples. Furthermore, the problem can 
be solved completely when the gauge dynamics is replaced by a Potts model. On the other hand in QCD with light 
quarks the solution will require a partial re-summation over both fermionic and gauge degrees of freedom. The 
new approach points to unexplored directions in the search for a solution to this more challenging sign problem. 
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1. INTRODUCTION 

Although lattice field theories allow for a com- 
pletely non-perturbative definition of a quantum 
field theory, computational limitations impose se- 
vere restrictions on our ability to perform calcu- 
lations in them. At the beginning of the new 
millennium many important questions that arise 
in a variety of fields like strongly interacting 
fermionic systems, frustrated systems, systems 
with (9-vacuum, chiral gauge theories and real 
time dynamics in quantum field theories, remain 
unanswered. A common problem that hinders 
progress in all these areas is the "sign problem" , 
which arises due to the inability to isolate the ap- 
propriate probabilistic measure for a stochastic 
computation. 

One is typically interested in the partition func- 
tion of a quantum many body system which can 
be compactly written as 



Z = Tr[cxp(-tf/T)] = £V[<7] 

[C] 



(1) 



where the sum over configurations C is used to 
evaluate the trace through the Feynman path in- 
tegral in a convenient basis. In general W[C], 
the Boltzmann weight of the configuration, is not 
guaranteed to be positive definite. The difficulty 
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in finding the right class of configurations such 
that W(C) is positive definite is referred to as 
the sign problem. However, since Z itself is posi- 
tive definite a summation over a sufficiently large 
number of configurations should always yield a 
positive answer. Unfortunately, the number of 
configurations required grows exponentially with 
the volume; it is hopeless to require even a com- 
puter to perform the re-summation. 

Recently, a lot of progress has been made in 
solving sign problems using cluster methods since 
it is sometimes possible to analytically sum over 
all the configurations that are obtained through 
cluster flips. Furthermore, this class of configu- 
rations can grow exponentially with volume. In- 
terestingly, every known type of cluster algorithm 
has been applied to solve to a new type of sign 
problem. The first solution was found in a classi- 
cal 0(3) model with a 9 parameter which makes 
the action complex JO . It was shown that a sum- 
mation over flips of Wolff clusters || produced 
only positive (or zero) weights. More recently, it 
was discovered that in certain four-Fermi mod- 
els a similar result emerges when one uses loop 
clusters known from quantum spin systems ||. 
This has been used to eliminate sign problems 
that arise due to Fermi statistics . Sign prob- 
lems in Potts models due to a complex action can 
be solved using clusters of the Swendsen-Wang 
algorithm ||. 

In this article the above progress is reviewed us- 
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ing the example of the Potts model which relevant 
to QCD at a finite density of static quarks. The 
next section contains a review of the sign problem 
in QCD with a new perspective on its "physical 
origin". It is argued that the sign problem has 
two parts, one coming from Gauss law type con- 
straints and the other from fermion permutations. 
The later sections show that the recent progress 
has produced examples in which both types of 
sign problems can be solved separately. In simple 
cases it may be possible to find solutions to the 
combined problem in the near future. 

2. THE SIGN PROBLEM IN QCD 

The partition function in QCD is usually writ- 
ten as a path integral over fermionic and gauge 
fields. Instead of representing fermions as Grass- 
mann variables, for computational reasons they 
are integrated out in favor of a fermion determi- 
nant. The final result takes the form 

Zqcb = J [dU] exp(-5 G [C/]) Det(D[U}) (2) 

where D is the Dirac operator. Given that the 
gauge field action Sg[U] is real the usual claim 
is that sign problems arise from the fermion de- 
terminant. As is well known in the presence 
of a chemical potential the fermion determinant 
Det (D[U]) becomes complex. Although this de- 
scription is quite accurate, it is a bit unsatis- 
factory since a lot of interesting dynamics of 
quarks has been pushed into the calculation of the 
fermion determinant. Exposing the anatomy of 
the fermion determinant may reveal some deeper 
"physical origin" to the sign problem. This can 
be illustrated by a simple example. 

Consider compact QED instead of QCD. The 
partition function in that case is similar to (||) ex- 
cept that the integration is over U(l) gauge fields 
instead of SU(3). Again adding a chemical po- 
tential makes the fermion determinant complex 
and hence introduces a sign problem for the same 
reasons as in QCD. But in this case it is possible 
to understand the origin of the sign problem at 
a more deeper level. If the fermion determinant 
is expanded as a sum over electron and positron 
world line configurations, the configurations that 
dominate at non-zero chemical potentials would 



contain more electrons than positrons. However, 
due to Gauss's law it is impossible to have more 
electrons than positrons in a periodic box. Hence, 
the weight of such configurations must be zero 
once Gauss's law is imposed. Since it is the in- 
tegral over the gauge fields that enforces Gauss's 
law, it is natural to have both positive and nega- 
tive contributions in the path integral for a fixed 
gauge field configuration. The negative signs can 
help cancel the wrong configurations. Under- 
standing the "physical" origin of the sign prob- 
lem makes the solution clear. In the case of QED, 
one must throw away all fermionic configurations 
in which there are more electrons than positrons. 
This effectively means the chemical potential has 
no effect on the partition function which is a well 
known fact JtJ. 

Unfortunately, fermionic world line configura- 
tions have sign problems of their own that arise 
from the Pauli principle. When fermions travel in 
imaginary time their positions can get permuted. 
Whenever this permutation is odd the configura- 
tion weight is negative. Thus, performing the in- 
tegration over the gauge field configurations alone 
in general will not yield a positive weight for a 
fermionic configuration. In fact the sign prob- 
lem in QCD arises due to an intricate mixture 
of fermion permutation sign and a more difficult 
constraints one of which is the Gauss's law. The 
example from QED suggests that it may be use- 
ful to consider a regrouping of both fermion and 
gauge field configurations in the search to find 
positive definite weights. 

3. THE STATIC QUARK LIMIT 

There are special limits where the sign prob- 
lem becomes much simpler. One is the limit of 
static quarks. Here the fermions cannot permute 
with each other and so problems arise only from 
Gauss's law type constraints. This simplifica- 
tion has prompted a variety of research on the 
physics of QCD with a finite density of static 
quarks. However, adding chemical potential to 
quarks that are infinitely massive is a bit subtle. 
In order to induce a finite density of static quarks 
one needs to take the chemical potential to infin- 
ity along with the quark mass always keeping the 
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chemical potential larger than the rest mass of the 
quarks. This limit was originally formulated in 
and later studied in [||. More recently the par- 
tition function for a hxed baryon number in the 
static quark limit was studied in ||. An equiv- 
alent but simpler approach to the limit of static 
quarks can be formulated by recognizing that the 
partition function for n quarks and fi antiquarks 
of mass m in the static limit is given by 



[dU] 



-S G [U]-m(n+n)N t 



i.l 



(3) 



where SgI^] is the gauge action on a N^N t lat- 
tice. The gauge coupling (3 is implicit in the defi- 
nition of Sg- The quantity $ is the sum over the 
trace of Polyakov lines at every spatial point x, 



(4) 



with Uq(x, t) representing the time like links. The 
grand canonical partition function is then given 
by 



Z = Z n , n exp -n)/T) 



[dU] e 



-S G [U]+h <i>+h' <S>* 



(5) 



(0) 



where h = e -( m -^ Nt and h' = e -^ m+ ^ Nt , The 
static quark limit is obtained by taking the m — > 
oo and /i — > oo limits with h = e~( m ~^ Nt fixed. 
In this limit h' — ► 0. Note that although the 
above approach treats quarks as bosons instead 
of fermions, at finite densities in the continuum 
limit it is impossible to distinguish between them. 

Due to the complex nature of e' 1 * the parti- 
tion function of eq. (JsJ) suffers from a sign prob- 
lem albeit a much simpler one as compared to eq. 

. Absorbing the magnitude of this term in the 
Boltzmann weight, one can define the average of 
the phase as, 



(Phase) = 



f[dU] 



-S G [U] + hRc<S> i/ilm* 



(7) 



Since (Phase) is a ratio of two partition functions, 
in the large volume limit it will vanish exponen- 
tially. In order to make progress without solv- 
ing the sign problem it is necessary to remain in 
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Figure 1. The phase diagram of static QCD. The 
point C is expected to be in the universality class 
of the 3-d Ising model. 



volumes where the average sign is not too small. 
This usually puts severe restrictions on our abil- 
ity to address interesting questions regarding the 
phase transition. 

The phase diagram in the f3 — h plane was first 
studied in |8| . Figure |l] shows the most likely sce- 
nario based on the well known fact that there is a 
first order transition on the h = line at (3 C due to 
the spontaneous breaking of the Z3 center sym- 
metry. When h ^ the symmetry is explicitly 
broken. However, the first order phase transition 
should typically persist for non-zero h with de- 
creasing latent heat ending at the second order 
critical point C which should be governed by the 
3-d Ising universality class. Although, the study 
in H could not confirm this picture it showed that 
the point C must be quite close to the h = axis. 
In particular the well known discontinuity in the 
Polyakov line at f3 c disappeared at the smallest 
h studied. This information suggests that with 
out solving the sign problem it may still be pos- 
sible to study the behavior of the model close to 
the critical line of the phase diagram on lattices 
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whose size is constrained by the closeness of the 
point C to the h = axis. 

4. SOLVING THE SIGN PROBLEM 

An interesting challenge is to be able to solve 
the sign problem arising in eq. (Q). Such a solu- 
tion can perhaps teach us about solutions to other 
more interesting sign problems in gauge theories. 
In order to make progress it is useful to learn from 
a simpler example where the gauge dynamics is 
replaced with a simple 3-state Potts model. The 
new partition function takes the form of 

Z = V" e 13 E[z]+hM[z] + h' M' [z] (g) 

Z 

with the energy and magnetization defined by 

e[z) = Y, 6 *~w m w = ( 9 ) 

x,i x 

Due to its symmetry properties, the three state 
Potts model with z x G {1, e i2jr / 3 , e~ i27T/3 }, has 
been considered a useful effective model for static 
QCD close to the phase transition JTc| , [TT| ] . 

An interesting solution to the sign problem in 
(||) was proposed in 113). It was shown that the 
partition function can be rewritten in terms of 
new variables which describe classical statistical 
mechanics arising from the Hamiltonian, 

H = J2 a \ M + Z)( M I "xl - /i fix), (io) 

x,j X 

that represents the energy of configurations 
labeled by quark number distribution n x G 
{—3, —2, 3} associated with sites and color flux 
variables G { — 1,0,1} associated with links. 
The configurations further obey the Gauss's law 
constraint 

3 

^('x,,-ix-,)=Hxmod 3. (11) 

This type of flux tube models is well known in the 
literature Jl3| . Since the Hamiltonian is real, ex- 
cept for the complications arising from the Gauss 
law constraint, the sign problem is completely 
solved. The introduction of new variables has 



helped in re-summing over classes of configura- 
tions in the original variables. In jl2| a metropolis 
algorithm was used to study the flux tube model. 

There is another way to solve the sign problem 
arising in (|J) based on cluster algorithms for Potts 
models ||. The essential idea in this case is to 
write the weight associated with a bond, 

e «s°*,*v = J2 [^,i^^(e K -l)+<5 M ] , (12) 

6=0,1 

as a sum over new bond variables 6, with 6=1 
representing the presence of a connection between 
the sites and 6 = representing the absence of a 
connection. Every configuration [z, b] of bonds 
and Potts spins is made up of a connected clus- 
ter (C) of sites, with the property that each con- 
nected cluster of sites always contains the same 
Potts spin zq. This property allows the Boltz- 
mann weight of each configuration [z, 6] to be 
written as a product over cluster weights and the 
partition function takes the form 

Z = E 1 1T W ^ C ) ex P ( L c [hzc + tiz* c \) 1 ,(13) 

[2,6] I C J 

where Lc represents the size of the cluster C and 
W(C) = (e K - l) Bc where Be stands for the num- 
ber of 6 = 1 bonds in the cluster. Since W(C) 
does not depend on the spin variables, it is easy 
to perform a sum over the allowed cluster spin 
variables z c G {l,e i2jr / 3 and e~ i27r / 3 }. The parti- 
tion function can finally be written as 

Z = E jll^) ^e^^+^j j. (14) 

It is easy to check that the Boltzmann weight in 
( |l4| ) for h,h' > is always positive which solves 
the sign problem. 

The two solutions to the sign problem sketched 
above suggest that an appropriate regrouping of 
the gauge configurations may yield a solution in 
the case of QCD with static quarks. However, 
this has not yet been achieved and remains an 
open problem for the future. On the other hand, 
based on the strong coupling expansion, a solu- 
tion of the type found in the flux representation 
of the Potts model may be easy to find. A cluster 
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Figure 2. Expectation value of Re(M) as a func- 
tion of (3 for h = 0.0005 (open circles), 0.001 (di- 
amonds), 0.0015 (filled circles), 0.002 (squares), 
and 0.0025 (triangles) in the infinite volume limit. 
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Figure 3. The critical line in the (3 — h plane. The 
data points represent the infinite volume extrapo- 
lations of the peaks of the specific heat. The solid 
line is a plot of eq. (H). 



type solution, although more exciting, could be 
difficult. 

5. RESULTS IN THE POTTS MODEL 

Recently, the partition function of eq. (||) has 
been studied numerically for h — h! . Since in 
this case the action is real the model does not 
suffer from a sign problem. The critical region 
has been analyzed in great detail and the critical 
end point C is found to be (J3 C = 0.54938(2), h c = 
0.000516(7)) QJll- The reader is referred to the 
original article for further results and details. 

The model with h' — is more closer in spirit to 
the physics of a finite density of static quarks. Al- 
though in this case the action is complex, the sign 
problem can be completely solved as discussed in 
the previous section. Using this solution and an 
algorithm described in |lq| , the behavior of the 
magnetization M was studied on lattices as large 
as 40 3 . Using finite size scaling analysis the were 
extrapolated to the infinite volume limit. In fig- 
ure H some results are plotted as a function of (3 
for various values of h. Clearly, even for h = 0.002 
the transition has already turned into a smooth 
cross over. This result is quite consistent with the 
results of ||. 



Since the critical line is very close to the h = 
axis, the sign problem is rather mild. An estimate 
shows that (Phase) ~ exp(— TV^/73 3 ) suggesting 
that only for N s > 80 will the solution to the sign 
problem become necessary. In fact it appears to 
be more efficient to simulate the theory with the 
real action and absorb the complex phase into 
observables for smaller lattice sizes. Using this 
technique simulations were extended to 80 3 lat- 
tices and the location of the critical line was de- 
termined using peaks in the specific heat. Figure 
|[ shows the data along with a fit to the critical 
line 

(3 = 0.550557(21) - /i/0.4337(78). (15) 

This line extrapolates very nicely to the first or- 
der critical point on the h — axis, which is well 
known to be f3 c = 0.550565(10) @. A prelim- 
inary analysis of the critical behavior near the 
end point C shows consistency with the 3-d Ising 
model. The critical point is found to be shifted 
mildly to (p c ~ 0.5495(1), h c ~ 0.000468(5)). A 
more detailed discussion of the analysis will be 
given in Jl6[ ] . 
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6. FERMION SIGN PROBLEMS 

The difficulties associated with finite density 
QCD are shared by much simpler four-Fermi 
models like the repulsive Hubbard model away 
from half filling. The conventional approach in 
interacting fermion models is to rewrite the inter- 
actions as fermion bilinears using auxiliary fields 
so that the fermions can be integrated out in favor 
of a fermion determinant, very much like in QCD. 
Hence it is not surprising that similar problems 
arise even in simple cases. An interesting question 
is whether one can deal with fermion interactions 
differently. Perhaps a more direct approach in- 
volving fermionic configurations will yield novel 
solutions if the fermion sign problem is tackled. 

Some years ago this approach was investigated 
in jig}] . Starting with finite density QCD in the 
strong coupling limit the gauge integration was 
performed explicitly and a purely fermionic ac- 
tion was obtained. The fermionic theory had 
six-Fermi couplings representing a baryonic mass 
term. Although one would naively disregard this 
approach as difficult, the fermionic model could 
be solved with a novel algorithm which used a 
clever partial solution to the sign problem in the 
model. Unfortunately, the physics of the strong 
coupling limit appeared rather uninteresting at 
high densities. It would be interesting to inves- 
tigate this approach further and ask if one can 
extend the study beyond the leading strong cou- 
pling approximation. 

Over the last couple of years a new class of so- 
lutions to sign problems in a variety of four-Fermi 
models have emerg ed The models are for- 

mulated directly in the fermionic Fock space, and 
configurations are labeled by fermion world lines. 
When the fermion permutation sign is taken into 
account with other local negative signs arising 
from transfer matrix elements, the partition func- 
tion can be written as 

Z = Sign[n] W[n] (16) 

[n] 

where [n] represents a fermion occupation number 
configuration. Since the Hubert space of fermion 
occupation number states is isomorphic to that 
of a spin-1/2 particle, the cluster algorithms for 




1 2 3 4 5 

Figure 4. A typical configuration of fermion oc- 
cupation numbers and clusters. The filled and 
open circles represent n — 1,0. The dark lines 
represent connected bonds. Shaded regions rep- 
resent interactions. There are five clusters in the 
configuration. 



quantum spin systems can be used to formulate 
algorithms for the fermions fisfl . Except for sub- 
tleties arising due to Fermi statistics, this essen- 
tially means that one is able to rewrite the parti- 
tion function (|l6| ) in terms of fermion occupation 
numbers and bonds 6, 

Z=53Sign[n,6]W[»,6], (17) 

[n,5] 

which is analogous to the example in the Potts 
model discussed in section 4. Figure || shows a 
typical fermionic configuration with an odd per- 
mutation along with clusters which have the prop- 
erty that by changing all the fermion occupation 
numbers within a cluster from to 1 and vice- 
versa, one obtains another allowed configuration. 
This update is referred to as a cluster flip. 

Recently, a relation between the change in 
the fermion permutation sign due to a cluster 
flip and the cluster topology was discovered pTJ] . 
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The relation showed that in certain models the 
sum over Boltzmann weights of configurations ob- 
tained through cluster flips for a fixed set of clus- 
ters is always positive; a result strikingly sim- 
ilar to the example of the Potts model. This 
observation has been used to design an efficient 
fermion cluster algorithm and is referred to as the 
"meron cluster algorithm" . It has been applied 
extensively to the study of the critical behavior 
near a Z2 finite temperature chiral phase transi- 
tion with staggered fermions 21 2^). This study 
for the first time demonstrated the emergence of 
Ising critical behavior in a fermionic model con- 
firming dimensional reduction. Figure |^ shows 
the chiral condensate as a function of tempera- 
ture in the staggered fermion mode. The scaling 
form A(T C — 7 1 )°- 314 ( 7 ) consistently describes the 
data which is with in two standard deviations of 
the behavior in the 3-d Ising model |22| ]. Dimen- 
sional reduction in fermionic models have been 
questioned in the past based on large N calcula- 
tions and results from hybrid Monte Carlo algo- 
rithms J23}] . The studies using meron algorithms 
are beginning to convincingly show the validity of 
the conventional picture. 

Although the number of models that can be 
solved with a meron cluster algorithm is limited, 
it continues to grow with time. Recently, new 
models in the family of the attractive and re- 
pulsive Hubbard models were found to be solv- 
able. It is possible to add a chemical poten- 
tial in the attractive model without violating the 
properties necessary for the meron algorithm to 
work. Such models are solvable with the hybrid 
Monte Carlo algorithm. However, previous stud- 
ies were limited to small sizes [ p4j . The meron 
algorithm on the other hand can be implemented 
on much larger lattices with relative ease. Re- 
cently, this has allowed a high precision study of 
the Kosterlitz-Thouless transition. More details 
of this work can be found in |^5| . 

7. FUTURE DIRECTIONS 

The amount of progress over the last cou- 
ple of years has been remarkable. A variety of 
physics associated with the spontaneous break- 
ing of continuous and discrete symmetries in 
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Figure 5. The critical behavior of the chiral con- 
densate as a function of temperature in a four- 
Fermi model. The open squares represent results 
from a meron algorithm on a 32 3 lattice at a mass 
of m = 0.001. The black circles represent the infi- 
nite volume extrapolation. The complete results 
are presented in p3|. 



certain fermionic models can now be studied 
with high precision using a meron-cluster algo- 
rithm. Questions related to chiral perturbation 
theory, physics of resonances, universality classes 
of phase transitions, etc. which appeared difficult 
with hybrid Monte Carlo, can be studied much 
more easily. It is important to extract as much 
as possible from these developments. 

In the context of finite density physics, there 
are many strongly interacting fermionic models 
that are still intractable. The example of QCD 
is a fundamental one. However, there are other 
models like the repulsive Hubbard model away 
from half filling or effective field theories of many 
nucleoli systems, that are equally exciting to 
tackle. Although the solution to the problem in 
QCD still appears distant, the new progress sug- 
gests that solutions to the simpler theories may 
not be too far away. However, one should to keep 
in mind the difficulty of formulating cluster algo- 
rithms since the most interesting solutions to sign 
problems have emerged in the context of such al- 
gorithms. Limitations in cluster techniques may 
be correlated to our limitations in finding solu- 
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tions to new sign problems. 

Based on the progress made so far many new 
questions for the future emerge. For example: 

(1) Can one find a meron algorithm in a model 
whose ground state supports long range 
fermionic excitations? This can open up the 
possibility of studying a variety of interest- 
ing quantum phase transitions. 

(2) Is it possible to find an algorithm where one 
can update both fermionic and bosonic de- 
grees of freedom together in an interacting 
model? 

(3) Can one extend cluster algorithms to new 
types of spin and gauge models? This 
may be possible by understanding the con- 
straints imposed by the topology of the con- 
figuration space pf| . 

The answers to these questions have the potential 
to determine if indeed the progress of the past 
couple of years is just the tip of an iceberg. 
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